Feshbach molecule production in fermionic atomic gases 
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This paper examines the problem of molecule production in an atomic fermionic gas close to an 
s-wave Feshbach resonance by means of a magnetic field sweep through the resonance. The density 
of molecules at the end of the process is derived for narrow resonance and slow sweep. 



The problem of the molecule production in an exper- 
iment where a system of fermionic atoms is tuned from 
far above to far below the Feshbach resonance has re- 
cently been studied both experimentally and theoreti- 
cally. First of all, experiments with the systems of ul- 
tracold atomic gases close to Feshbach resonance [H, ^ 
explored this physics. A number of theoretical papers 
followed 0, II E 0, 0, Q. In the language of the two 
channel model [9|] , the problem can be formulated in the 
following straightforward way. Consider the Hamiltonian 
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fermions (atoms) of momentum p and spin cr, b]^, 6q are 
creation and annihilation operators of bosons (molecules) 
with momentum q, g is the coupling, V is the volume, 
and eo is the parameter, called detuning, which controls 
whether fermions or bosons are energetically favorable. 
The system defined by Eq. ^ has been extensively stud- 
ied in the literature, and it is well know that for large 
positive eg this system is a BCS superconductor, while 
when eo is large and negative, this system describes the 
Bose condensate of weakly interacting bosonic molecules. 

Now suppose eo is a function of time (representing a 
sweep through the Feshbach resonance), given by 



(2) 



where A is the rate of its change. Furthermore, suppose 
initially, at t ^ — oo, the system sits in its ground state 
which is represented by the Fermi sea with the Fermi 
momentum pf in the absence of the molecules 



By dimensional analysis, the density of the created 
molecules nf, has to be proportional to the initial density 
of the fermions, np — pp/iQir^) {np is the density of each 
of the species, spin up or spin down). It can also depend 
on the two dimensionless parameters in this problem, one 
being the Landau-Zener parameter [TO, El 
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which controls the rate of the sweep through the Feshbach 
resonance. The other is the parameter which controls the 
width of the resonance [l2, 13, 14 1, 
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are creation and annihilation operators of jjere 
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is the Fermi energy of the initial fermion distribution. 

Thus rif, = npf(T,"f), where / is some dimensionless 
function which needs to be determined. In this paper 
we analyze the regime where 7 ^ 1, In [l/(r7)] <C T ^ 
1/7, corresponding to a relatively slow sweep of a narrow 
resonance. Our main result is 
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In Ref. Pi the molecular density rib was calculated for 
r <C 1 perturbatively in powers of F, equvalent to study- 
ing a fast sweep. The answer was given by 
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where 9{x) is the usual step- function, 9{x) = 1 if a; > 
and 9{x) = otherwise. One is asked to find the density 
of molecules at the end of the process, or compute 
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(notice that it is independent of 7 within this order of 
perturbation theory) . This should be contrasted with the 
naive ap plic ation of the Landau-Zener formula for level 
crossing [lOl . ITU (which requires some generally unjusti- 
fied replacement of Eq. ([T]) by some sort of a two-level 
system) which would have given 
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fermions, to replace Eq. pTjl with 



Clearly already in the second order of perturbation the- 
ory, the Landau-Zener formula breaks down. Thus a re- 
maining interesting question would be to compute rih at 
large F, the task beyond the reach of the perturbation 
theory and accomplished in this paper. 

We start our derivation of Eq. ([8]) with replacing the 
Hamiltonian Eq. ([T]) by a reduced Hamiltonian where 
only the bosonic zero (momentum) mode is kept 
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(11) 

This approximation, quote common in the literature ded- 
icated to this problem @, [1] , is inspired by the fact that 
one expects a finite fraction of bosons to be Bose con- 
densed at i ^ cxD, at least at slow enough sweep A. Thus 
perhaps only the zero mode of the bosons should con- 
tribute in a substantial way. 

We are going to argue that going from Eq. ([T]) to 
Eq. (|lip can only be justified if 7 ^ 1 and F ^ 1. In- 
deed, it is well known from the studies of the equilibrium 
problem (where is time independent) that Eq. (|lip is 
a good approximation to Eq. ([1]) if 7 <C 1, or if the res- 



onance is narrow i l3l . Il4l | . Additionally, at large positive 



eo, where the system described by Eq. ([T]) is in the BCS 
regime, Eq. pT|) is a good approximation regardless of 
the value of 7. On the other hand, if 7 3> 1, then in our 
problem where eo changes in time from large positive to 
large negative values, a range of values of Cq must be tra- 
versed, called the unitary regime, where Eq. (jlip breaks 
down. Thus it is natural to conjecture that, if 7 3> 1, 
Eq. pT|) is not a good approximation to Eq. ([1]) not only 
in case of the time independent problem but also for the 
time dependent sweep Eq. ([2]). 

At the same time, if 7 <C 1, then if F ^ 1 (very 
fast sweep) a straightforward perturbative calculation for 
Eq. ([TT|) shows that the density of the produced bosons 
goes as rif, « T/V, or in other words, in the thermody- 
namic limit it vanishes. This is quite different from the 
analogous result for Eq. ([J), given by Eq. This oc- 
curs because the molecules in Eq. ([T]) formed under the 
conditions of a fast sweep are not Bose-condensed. Thus 
in this regime Eqs. ([T|) and PT|) describe completely dif- 
ferent physics. 

This leaves the parameter range 7 <C 1, F ^ 1 as the 
only regime where Eq. ([1]) and Eq. (jlip may be equiva- 
lent. And indeed, since under a very slow sweep one ex- 
pects that most produced molecules are Bose condensed, 
it is not unreasonable to assume that this is the case. 
From now on we adopt this point of view. 

Subsequently, additional common approximation is of- 
ten made by further neglecting the dispersion of the 
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We would like to argue that reducing the problem to 
Eq. (I12|) is problematic. Indeed, it is well known that a 
time independent version of the problem Eq. (|12p un- 
dergoes a quantum phase transition as a function of 
eo [3 [3 ; ^ feature absent in Eq. ^ which is widely 
believed to go through a crossover (termed the BCS- 
BEC crossover in the literature) [l3]. Thus this model 
describes quite a different physics (see however Ref. 
where it is argued that at large 7, the large g makes it 
possible to neglect both fermionic and bosonic disper- 
sions and immediately arrive at Eq. ()12p . in fact in the 
sector where it does not have a phase transition). 

In this paper, we present the solution to the problem 
defined by Eq. (P) with Eqs. ^ and ©, with 7 < 1 and 
at a slow sweep F 1, thus justifying replacing Eq. ([1]) 
with Eq. (fTTj) . 

Next, following the standard approach, we would like 
to replace the fermion operators by the Anderson pseu- 
dospin operators defined by 
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We observe that within the model Eq. (fTTj) one can safely 
replace the quantum operator 6 by a classical field b. 
Ref. which worked with the Hamiltonian Eq. (|12p. 
examined this question in a lot of detail. It was estab- 
lished that the classical dynamics was indeed a good ap- 
proximation to the quantum behavior of Eq. (|12p with 
the exception of the early evolution of the system where 
the number of bosons was close to zero. A crucial dif- 
ference, however, between Eq. and Eq. (|lip studied 
here is the existence of a quantum phase transition as eg 
in Eq. changes. As a result, a pure classical variable 
b{t), if it is zero initially, remains zero forever within the 
classical dynamics of Eq. (|12p . and its quantum behavior 
must be studied to find nt = \b\ in the infinite future. 
Fortunately, this does not happen in Eq. (fTTj) . 

Thus let us replace the operator b by the classical vari- 
able b. Once this is done, the spin dynamics can be 
described by the classical canonically conjugate variables 



and do where 
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Here rip gives the occupation number of a fermion pair 
at momentum p. Taking into account that the Hamilto- 
nian Eq. (jlip conserves the total number of particles, or 
that X]q"-q + l^l = npV, wc find the effective classical 
Hamiltonian equivalent to Eq. pT|) 
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\/np{l - Up) cos (j)p, 



where a convenient notation is introduced, Cp = / (2m). 
This Hamiltonian is equivalent to the classical version of 
Eq. (jlip . and is supplemented by the initial conditions, 
set at t = — oo, stating that Up = 6{pf — p)- We could 
now take advantage of the integrability of Eq. (flS)) when 
Eq is time independent, discovered recently in Ref. [lit . 
The integrability allows in principle to construct its adi- 
abatic invariants. Then once eo depends on time, one 
could calculate how the adiabatic invariants change with 
time following [l^. This technique, although promis- 
ing, is technically involved and will not be pursued here. 
Instead, we will follow a related yet somewhat different 
approach. 

The structure of the motion described by this Hamil- 
tonian is elucidated by its fixed points, or points where 
dH/dup = dH/d<f)p ~ 0. These are given by 



1 



1 - 



Here we employ standard notations 



(16) 



(17) 



where A and /i have to be determined by solving the two 
coupled equations 



np 
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This is the standard BCS-BEC solution of the two chan- 
nel model, with Eq. (|18p representing the particle number 
and the gap equations. In fact, in the limit 7 — > 0, it is 
easy to solve these two equations 13|, |lj| , with the result 



eo > 2eF : A' = ^f, A = 0, 
2eF > eo > : ^ = eo/2, A = yCp" 
> eo : /i = eo/2, A = 



(2m^)3/2 



67r2 



/np- 



(19) 



Notice that if 7 > 0, then the solution deviates slightly 
from Eq. In particular, at large eo, A is no longer 

zero. Now if eo depends on time according to Eq. ([2]) with 
very small A, or F 1, then the mechanical system will 
adiabatically follow the stationary solution Eq. (fT6|) with 
Eq. (fT9|) . This corresponds to the complete conversion of 
the fermions into bosonic molecules, since the density of 
molecules at the end of the process will be given by A^ 
which will be equal, at large negative eo, to np. 

However, as the rate A is increased, the number of the 
produced molecules will decrease. We can describe this 
in the following way. Let us compare the Hamiltonian 



Eq. ([15)) with the standard Hamiltonian of the Landau- 
Zener two level crossing. That problem can be set up in 
terms of the two time evolution equations 



^ -\tiJi+ gtp2, = Xtilj2+ gtpi (20) 

with the initial conditions set at t —00, ipi — and 
-02 = 1- The Hamiltonian for this problem is given by 



.dtp2 
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LZ 



At(|V2|'-|V'i|') -)-g(V'i>2 + 02>i)- (21) 



Under the mapping ipi = \/l — ne"^S "02 = \/ne^'^^, 
<j} = <j)i — 4>2 this Hamiltonian can be transformed into 



LZ 



2Xtn + 2g\/ n(l — n) cos (j) 



(22) 



with the initial condition n = 1 at i = — 00, which pos- 
sesses a close similarity to Eq. ([T5|). The quantity of 
interest is n computed ai t ~ 00, known to be equal to 



At=oo = exp(-7r5VA) • 



(23) 



We can then interpret Eq. (|T5|) . on the assumption 
of 7 <C 1, in the following way. Each Up undergoes its 
own Landau-Zener transition when eo is close to 2ep with 



an effective coupling given by g^np — X)q"q/^- While 
this process goes through, np for other values of p do 
not change appreciably (only one value of p is in 'reso- 
nance' at a given time), and neither does riq. Within 
this interpretation, the final value of rip at large positive 
times, which we denote n^, can be found by comparing 
Eq. (HID with Eq. ^ and using Eq. 



(24) 



Here nh{p) has a meaning of the boson density at the 
point in time where eo = 2ep, and is given by 



(25) 



p<q<PF 



where riob is the boson density at eo = 2eF- 

Introducing a variable x — j (Ott^), which varies from 
to rii?, its easy to see that Eqs. (p4)) and (|25p are equiv- 
alent to the following equation 



dni,{x) 
dx 



1 



(26) 



which leads to the following boson density at the end of 
the process or when a; = 



Ub = 7if,(0) = In 



- 1 



1 



(27) 



Notice that this argument assumes ng'^nob/X <C 1. Oth- 
erwise, Eq. (|27p can predict that rib > np which is con- 
tradictory. The origin of this restriction lies in the fact 
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that at those very slow rates where this condition is vio- 
lated, one needs to take into account that the state of the 
fermions at eo = 2eF is a paired superfluid and not a filled 
Fermi sea. This leads to a rather complicated problem 
since Eq. ((24)) no holder holds true as initially were not 
1. Considering this goes beyond the approach developed 
in this paper. 

In the limit 7 — > 0, not is equal to zero. A direct sub- 
stitution of that into Eq. (P7)l leads to Ub = 0. That 
means, if we start with zero bosons, we end with zero 
bosons. This situation is reminiscent of what happens 
with Eq. (fT2| as discussed at length in Ref. [l6|. How- 
ever, as was pointed out earlier in this paper, our case is 
crucially different in that, unlike Eq. (fT2|) or even Eq. (fTTI) 
in the limit 7 = 0, at any finite 7 our system undergoes 
a crossover and not a phase transition. To be more pre- 
cise, during the initial time interval where eo > 2ep, the 
system will adiabatically follow the solution of Eq. (fT8|) . 
This evolution is very slow at 7 <C 1, and the deviation 
from adiabaticity during this time is negligible. The value 



of A'' at eo = 2e_F should be used for no;, in Eq. p7|) . 
at eo > 2eF is small but nonzero, and one needs to go be- 
yond the approximation involved in Eq. (|19|) to find what 
it is. That value can be calculated by carefully solving 
Eq. (UH]) at finite yet small 7. Following, for example, 
the techniques described in Ref. (in particular their 
equations (6.35) and (6.36)), we find that at eg = 2eF, 
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(28) 



where k = {ep — fi) /ep. Solving this equation with log- 
arithmic accuracy at small 7, we find 



nob = A^ 



3(2m)V 
327r2VeF 



np ~ np^. 



(29) 



We now substitute this value into Eq. under the fol- 
lowing conditions r7 <C 1 (which follows from the con- 
dition TT^^rxQ^/A < 1) and F > In [l/(r7)], and find the 
final answer Eq. ([8]) which completes the derivation. 

We see that the answer, first of all, implies that it is 
not as easy to reach the adiabatic behavior by reducing 
the rate A as it would have been for a Landau Zener level 
crossing, where the approach to the adiabatic behavior is 
exponential. Here the number of remaining excitations 
is linear in A, not exponentially suppressed, similarly to 
the examples considered in Ref. 20] . To covert an appre- 
ciable fraction of atoms into molecules, the parameter F 
must be at least bigger than ln(l/7), the value which be- 
comes larger (thus the rate smaller) as 7 decreases. We 
also notice that Eq. ([8j looks similar to some of the ex- 
pressions found in Ref. [51 , yet the criteria of applicability 
found in this paper are quite different. 

We also see that the answer Eq. (|S]) bears some similar- 
ity to the behavior of the reduced model Eq. (fT^ studied 



in Ref. The similarity appears to have its origin in 
the fact that Eq. p2|) undergoes a phase transition as eo 
is varied, while in the limit 7^0, Eq. pT|) also has a 
phase transition, as should be clear from Eq. Yet 
at finite 7, Eq. pT|) undergoes but a crossover, thus the 
similarity is somewhat limited. 

Finally, one is tempted to ask what this implies for 
the molecule production at large 7^1, equivalently 
for broad resonance, the condition under which most of 
the experiments are done and the majority of theoretical 
papers dealing with the subject are written. Here one has 
to admit that the methods developed in this paper are 
not designed to deal with this problem. We also remark 
that the well known and experimentally relevant papers 
d, 01 study a regime quite different from the one studied 
here (corresponding to €0 = 1/0 — 9{t)Xt) so that a direct 
comparison with their results is not possible. 
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